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A method for studying the causal structure of space-time 
evolution systems is presented. This method, based on a gen- 
eralization of the well known Riemann problem, provides in- 
trinsic results which can be interpreted from the geometrical 
point of view. A one-parameter family of hyperbolic evolution 
systems is presented and the physical relevance of their char- 
acteristic speeds and eigenfields is discussed. The two degrees 
of freedom corresponding to gravitational radiation are iden- 
tified in an intrinsic way, independent of the space coordinate 
system. A covariant interpretation of these degrees of free- 
dom is provided in terms of the geometry of the wave fronts. 
The requirement of a consistent geometrical interpretation of 
the gravitational radiation degrees of freedom is used to solve 
the ordering ambiguity that arises when obtaining first or- 
der evolution systems from the second order field equations. 
This achievement provides a benchmark which can be used 
to check both the existing and future first order hyperbolic 
formalisms for Numerical Relativity. 



I. INTRODUCTION 

The structure of Einstein field equations has deserved 
great interest since the very beginning of General Rel- 
ativity. It was early noticed that, by rearranging the 
order or partial derivatives, the principal part of the four 
dimensional Ricci tensor could be written as a sort of 
generalized wave equation 



2 R„ v = -□<?, 



d v T, 



(1) 



where the box stands for the d'Alembert operator acting 
on functions and we have written for short 



r M = gP 



T" pa = -Ox" 



(2) 



This opened the way to the use of harmonic spacetime 
coordinates (Ox" - = 0) in order to obtain an hyperbolic 
evolution system ]|-^| . 

By the middle of the past century, however, the inter- 
est focused in the relativistic Cauchy problem. The 3+1 
decomposition of the line element: 



ds 2 



-a 2 dt 2 

7y (dx i + (3 l dt) (dx j + f3 j dt) i, j 



1,2,3 (3) 



allowed one to express six of the ten second order original 
equations as a system of evolution equations for the met- 



ric 



and the extrinsic curvature Ka of the t = constant 



slices, namely: 



(dt 
(d t - 



-2 a K k 



(4) 

[^R ij -2K%+trK K i5 ] (5) 



where C stands for the Lie derivative (we restrict our- 
selves to the vacuum case for simplicity). The remaining 
four equations could instead be expressed as constraints: 



(3) R - tr(K 2 ) + (trKf = 
V k K k i-di{trK)=0 



(6) 
(7) 



This opened the door to a new way of obtaining hyper- 
bolic evolution systems The key point was to use 
in one or another way the momentum constraint (|7|) to 
ensure hyperbolicity while keeping the freedom of choos- 
ing arbitrary space coordinates on every t = constant 
slice. This allows for instance to use normal coordinates 
= 0) without affecting the mathematical structure of 
the evolution system. (See Ref. JhJ for a detailed compar- 
ison of the "old" and "new" hyperbolic systems). 

These findings came at the right moment for people 
working on Numerical Relativity with a view on the grav- 
itational waves detector projects starting by the turn of 
the century. Following the wake of these first works, 
many groups found their own way of combining the mo- 
mentum constraint with the evolution equation (^|), lead- 
ing in each case to a new brand of hyperbolic systems 

Suddenly, the problem of having too few hyperbolic 
formalisms for Numerical Relativity turned into the op- 
posite problem of having too many of them. Some of 
the works considered even multiparameter families of hy- 
perbolic systems 1 14 1^J2^] . Faced with the problem of 
choice, we think that the right question at this point is: 
what are hyperbolic systems for? 

There are many answers, of course, but we can get a 
hint from the physical point of view when we realize that 
all but one of the resulting characteristic speeds can be al- 
tered when playing with these arbitrary parameters (see 
the work of the Cornell group for a fairly complete 
analysis) . The only one that has an intrinsic value (light 
speed) is associated with the two (transverse traceless) 
degrees of freedom of gravitational waves. The evident 
question now is wether or not we can directly relate the 
corresponding explicit eigenfields with gravitational ra- 
diation propagation. We will see below that a positive 
answer actually determines one of the parameters, intro- 
duced in the Cornell paper [^6). We claim that the re- 
quirement of having explicit eigenvectors describing the 
gravitational waves degrees of freedom in an intrinsic way 
(for an arbitrary choice of space coordinates) is a strong 
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benchmark to discriminate between different hyperbolic 
evolution systems. 

In order to obtain our results, it has been crucial the 
use of an intrinsic method in order to find the charac- 
teristic speeds and the corresponding eigenfields. This 
method, based on a General Relativistic analogue of the 
well known Riemann problem of Computational Fluid 
Dynamics (CFD), is completely independent of the space 
coordinate system used at any given t = constant hyper- 
surface. It is presented in Section 3 and then applied in 
Section 4 to a new simple one-parameter family of hy- 
perbolic systems in order to illustrate our arguments. 



II. A FIRST ORDER EVOLUTION SYSTEM 



Note however that the actual form of the fluxes asso- 
ciated to Kij in eq. (||), will depend on the expression of 
the three-dimensional Ricci tensor ( 3 'Hij m terms of the 
metric derivatives (see for instance Ref. [^6| for a discus- 
sion of this point). We will take in this paper the original 
expression, namely: 

^Rij — dkT k ij — djT k ki + r fe fe r T r ij — T h ri T r k j (15) 

This is in contrast with the actual choice in most of 
the hyperbolic formalisms published up to now, where 
a three-dimensional version of the decomposition (^) is 
used instead of (15) in order to mimic the properties of 



generalized wave equations. 



The evolution system (|4||5[) is of first order in time but 
second order in space. To obtain a system which is also of 
first order in space, we will follow the standard procedure 
by considering the first space derivatives of the metric 
coefficients 7^ as new independent variables: 



Dkij = 1/2 dkjij 



(8) 



The evolution equation of these new variables D ki j can 
then be taken to be: 



d t Dkij + d k (a K i: 







(9) 



where we are using for simplicity normal coordinates, so 
that (3 l = 0. The relationship (§p can then be understood 
as a first integral of the extended system (Q, ^, ^). 

In order to complete this system, however, we need to 
provide evolution equations for the lapse function a and 
its first derivatives. It can be done in a straightforward 
way as follows: 



dt In a = —a Q 
Aj = di In a 
d t Ai+di(a Q) = 



(10) 

(11) 
(12) 



where, again, the equation ([lT]) must be considered as 
a first integral of (|l2|). Notice that one could instead 
prescribe a more general evolution equation for Af. 

dtAi + d^a Q) + n a (V fc K\ - d, trK) = (13) 

where fi is an arbitrary parameter. In that case, ( |ll| ) 
would be a first integral of ([l3|) if and only if the momen- 
tum constraint (^) is satisfied. 

The complete first order evolution system can then be 
expressed in a balance law form, namely: 



d t u + d k F k 



(14) 



where u is the array of independent variables (the met- 
ric coefficients and their first derivatives), and both the 
Fluxes F k vector array and the Sources S scalar array 
depend of u. 



III. THE GENERAL RELATIVISTIC RIEMANN 
PROBLEM 

The balance law form ( |l4|) of the evolution system is 
familiar from the Computational Fluid Dynamics (CFD) 
domain, where we can consider (14) in the sense of distri- 
butions, so that discontinuous solutions for u are allowed 
(weak solutions). The discontinuity surface £ can be 
given by: 

<f>(x l ,t) = constant 
d t <j> + v* drf = , 

where v l is the velocity of the x l — constant points of 
the surface S. The consistency of these weak solutions 
with (|l4|) requires the exact cancellation of the Dirac 
8 terms arising where the derivatives of discontinuous 
functions are performed. One gets then the well known 
Rankine-Hugoniot conditions: 



v[u] =n k [F k ] 



(16) 



where the square bracket stands for the jump of the quan- 
tities across E, n k is the unit normal to £ and v = v l n% 
is the propagation speed of the discontinuity front. The 
eigenvalue problem ( jig ) is usually known in CFD as the 
"Riemann problem" . 

In General Relativity, we are familiar with 'matched' 
metrics having discontinuous first order derivatives and 
other less common weak solutions, like colliding waves. 
A peculiar feature in General Relativity is however that 
the non-linear terms are limited to the Sources S in p4|) , 
which can not originate 6 terms and do not appear then 
111 |t|). This means that the General Relativistic Rie- 
mann problem is a linear one and can be solved explic- 
itly. This implies also that the eigenvalues v of (|l|) will 
coincide with the characteristic speeds, contrary to what 
happens in CFD, where one finds a much richer structure 
so that, in addition to plain discontinuities, one can get 
also shocks and rarefaction waves propagating at non- 
characteristic speeds. We will take advantage of these 
simplifications in what follows. 
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Let us then write down the Riemann problem corre- 
sponding to our evolution system ||, [| |l|), taking 
from granted that the metric coefficients are continuous 
so that the only jumps can be on the first derivatives 
quantities, namely 



v[Kij] = a n k [\ k ij] 
v[Dkij] = a n k [Kij] 
v[A t ] = an k [Q S k i 



+ f i(K k l -trK 6 k i)} 



(17) 
(18) 
(19) 



where we have written, allowing for (|l5|), 

A fe y = -r\- + 1/2 S k (Aj + Dj) + 1/2 S k (Ai + A) (20) 

and D k = 7 y D kij . 

Let us remark that the Riemann problem ( |l7| - |l9| ) is 
stated for a generic orientation of the unit normal n k , in 
a way that is independent of the space coordinate system. 
In this sense, it provides an intrinsic method to obtain the 
characteristic speeds and the corresponding eigenhelds in 
terms of geometrical objects related with both the t — 
constant hypersurface and the characteristic surface X. 



GETTING AN STRONGLY HYPERBOLIC 
SYSTEM 



IV. 



We will now solve the eigenvalue problem (|17j-^9|) in 
order to display the causal structure of the evolution sys- 
tem (E, 0, H |l3|). We remember that the system will be 
strongly hyperbolic if and only if all the eigenvalues v are 
real and if the set of eigenhelds is complete, in the sense 
that it spans all our variable space. It follows easily from 
(U) that 



v[D ±ij ] = 



(21) 



where the symbol _L replacing one index means that we 
are taking the corresponding components tangent to X 
(orthogonal to n k ). The remaining degrees of freedom in 
( |l8| ) can then be written as 



v[D n ij\ = a[Kij] 
where we have noted for short D n i 



(22) 



n k D k 



Let us now analyze the "mixed" components of ( |l7| ) 



v[K n ±}= a[-D± nn 
v[A±] = a u[K n x ] 



1/2 A A 



(23) 
(24) 



so that, allowing for (|2lj), one gets the following set of 
eigenvectors 



K n x ± (1/2 A± — D x nn ) 



(25) 



with characteristic speeds v = ±a yj /Lt/2, respectively. 
This means that our evolution system will be (strongly) 



hyperbolic only if n > 0. This surprising parameter- 
dependence of the characteristic speeds can be under- 
stood if we recall that physical solutions must verify the 
momentum constraint (Q). For a discontinuous solution 
one gets 



[K n i\ = m[trK] 



(26) 



so that the first term of (|25|) comes from the transverse 
part of (p6|). This means that the characteristic cones 
spanned by (^5|) are just mathematical artifacts in or- 
der to get an hyperbolic system: they can not propa- 
gate physical information. We will take n = 2 to impose 
this arbitrary characteristic speed to coincide with light 
speed. 

Let us focus now on the purely transverse components 
in @. We get 



v[K xx ] = a[-T n 1 
v[D n xx ] = a[K xx ] 



so that, allowing again for (]2 
eigenfields 



J = a[D'\ ± - 2 D 



(-L-L) 



(27) 
(28) 



one gets the following 



K xx ± (D r 



2D 



(-L-L) 



(29) 



with characteristic speed v — ±a, independently of the 
choice of fi. Notice that the longitudinal component of 
(Eq) can be written as 



( 7 «_„< n j )[K ij ]=0 



(30) 



so that, again, the trace component in (|29|) can not de- 
scribe the propagation of any physical information. See 
for instance Ref. |2(J to verify that the corresponding 
characteristic speed can be also modified at will by using 
the energy constraint in a suitable way. We can conclude 
that only the traceless part of the transverse components 
( p9| ) can actually propagate physical information. These 
two degrees of freedom are then good candidates to de- 
scribe intrinsic features, like gravitational waves, as we 
will discuss further in the next section. 

The remaining degrees of freedom in the original sys- 
tem ( |l7Hl^ ) are given by: 



v[trK] 
v[A n ] 



a[A n + 2 (D n - D k 
a[2 (K nn - trK) + Q] 
a[trK] 



(31) 
(32) 
(33) 



where now we must specify the gauge-related function 
Q in terms of other quantities. We have postponed the 
gauge specification up to that point to emphasize that 
our previous analysis is gauge independent. Now we will 
take the simple prescription 

Q = f trK (34) 

where / can be again an arbitrary parameter. We get 
then, allowing once more for (|2l|): 
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v[A n + 2 (D n - D k kn ) -fD n }=0 



(35) 



and also that the following gauge-dependent combina- 
tions 



y/f trK ± (A n + 2 D" - 2 D k kn ) 



(36) 



are eigenfields with characteristic speeds v = ±a y 7 /, so 
that the system can be strongly hyperbolic only if / > 0. 
We can take for instance / = 1 if we insist in having light 
speed as a characteristic speed also here. 

The analysis is complete now, because for /i = 2, / > 0, 
the eigenfields @, (§§), @, @ and (§§) span the 
whole space of variables. It follows that the resulting 
evolution system is strongly hyperbolic, no matter of the 
particular orientation n k of the characteristic surface E. 



V. GEOMETRICAL AND PHYSICAL 
INTERPRETATION 

Let us now decompose the line element of the three- 
dimensional t = constant hypersurfaces in the following 
way (2+1 decomposition): 

7ij dx i dx 3 = -N 2 dz 2 

+ a ab {dx a + A Q dz) (dx b + X b dz) a, b = 1, 2 (37) 

where the coordinates x a span the characteristic surface 
E, cr ab = lab is the induced metric on E and z is a trans- 
verse coordinate so that the unit normal to E can be 
expressed as 

m = NSi* . (38) 

The eigenfields ( ^9|) can then be written as 

K ab ±N (D z ab -D ab z -D ba z ) (39) 

and a straightforward calculation shows that the extrinsic 
curvature of the two-dimensional surface E, namely 



n ab = 1/(2N) [d z a ab - C x (a ab )} 



is given by 



n ab = N (D z ab - D ab z - D ba z ) 



so that the eigenfields (29) can be expressed as 
K ab ± K ab 



(40) 
(41) 
(42) 



which shows that they are two-dimensional tensors as- 
sociated to the characteristic surface E in an intrinsic, 
coordinate-independent way. 

Looking for a physical interpretations, let us note that 
the traceless part of (^2|) contains the shear degrees of 
freedom as measured by a two-dimensional array of ob- 
servers sitting on the surface E. This is precisely the 
kind of effect one expects from gravitational waves with 



wavefront E. Our results can then be interpreted as an 
extension of the standard description of the effect of a 
gravitational wave on an array of freely falling (geodesic, 
a = constant) observers. This is an extension in the 
sense that our gauge conditions (Q) apply to more gen- 
eral kinds of observers (like static ones, where both Q and 
trK vanish separately so that (|34|) holds for any value of 
/). This can be relevant for modelling Earth-based grav- 
itational wave detectors, which are certainly not in free 
fall. 

Notice, however, that the eigenfields expression (^9|) is 
sensitive to the ordering ambiguity of space derivatives 
when passing from (||) to a first order system (see Ref. 
p6[ for details). The tensor character of the eigenfields 
is completely lost if one starts with a three-dimensional 
version of (^) instead of (^|). The intrinsic information 
resides of course in the system ( [Til) , where the divergence 
of the fluxes is taken, but the local covariance of the 
fluxes themselves is lost if we make a choice different 
of (|l5|). We claim that our results provide a criterion for 
solving the ordering ambiguity: the eigenfields describing 
the gravitational radiation degrees must have an intrinsic 
geometrical meaning independent of the space coordinate 
system. This opens the way to a coordinate-independent 
local description of gravitational wave detectors at Earth- 
based laboratories. 



VI. CONCLUSIONS AND OUTLOOK 

We have started this paper by proposing a new sim- 
ple method, based on the General Relativistic Riemann 
problem, in order to analyze the causal structure of any 
first order flux-conservative evolution system in an in- 
trinsic way, independent of the space coordinates. To 
compare with previous results jl2), we can look at the 
characteristic speeds and see how we get now light speed 
as v = ±a (length per unit coordinate time), whereas 
in coordinate-dependent analysis, the speed along the z 
axis is given by v = ±a ^/j zz (coordinate displacement 
per unit coordinate time). 

The use of this intrinsic method has allowed our- 
selves to explain the appearance of arbitrary parameter- 
dependent characteristic speeds p(| as a direct conse- 
quence of the mixing of evolution and constraint equa- 
tions. As we emphasized in Section IV, only two degrees 
of freedom are intrinsically related with causal propa- 
gation of physical information. Hyperbolic formalisms, 
which can only consist of characteristic cones and lines, 
apply instead to an extended, unconstrained, space of so- 
lutions. The subspace of physical solutions is recovered 
by restricting the initial conditions to those that satisfy 
the constraints. The extra degrees of freedom are then to 
be regarded as mathematical artifacts devised to ensure 
hyperbolicity and their corresponding (arbitrary) char- 
acteristic speeds are irrelevant from the physical point of 
view. 
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Regarding the two remaining degrees of freedom, that 
can consistently be interpreted as describing gravita- 
tional radiation, we have provided an intrinsic geomet- 
rical interpretation of the corresponding eigenfields in 
terms of the shear of the congruence of eulerian (labora- 
tory) observers combined with the shear of the front- wave 
surfaces S. The requirement that the eigenfields admit 
this geometrical interpretation allows one to solve the or- 
dering ambiguity for space derivatives when passing from 
the second order form of Einstein evolution equations to 
a first order system. 

The fulfillment of this requirement must be regarded as 
an important benchmark to check any hyperbolic formal- 
ism. This is actually a strong requirement: none of the 
hyperbolic formalisms proposed previously by the Palma 
group |||l^,|l3| verifies it. It is not difficult, however, to 
modify the existing formalisms to comply with it |p7fl . 

This "explicit gravitational radiation" requirement is 
mandatory if we want to provide a local description of 
gravitational radiation in order to model for instance the 
process of detection of these waves by a given array of 
observers. In this sense, our results generalize the well- 
known textbook results for freely falling observers, based 
in the geodesic deviation equations. In particular, our 
results can be useful to model Earth-based laboratories, 
where the Earth gravitational field can be approximated 
by the well known Schwarzchild solution and the other 
effects, including the gravitational wave itself, can be in- 
cluded using perturbation theory. 
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